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Abstract. We explain why the conventional argument for deriving the time-dependent Born-Oppen- 
heimer approximation is incomplete and review recent mathematical results, which clarify the situation 
and at the same time provide a systematic scheme for higher order corrections. We also present a new 
elementary derivation of the correct second-order time-dependent Born-Oppenheimer approximation 
and discuss as applications the dynamics near a conical intersection of potential surfaces and reactive 
scattering. 
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1. Introduction 

m 

■r^ • Through the discovery of Schrodinger the theoretical physics community attained a powerful tool for com- 
, puting atomic spectra, either exactly or in perturbation expansion. Born and Oppenheimer [4] immediately 

■ strived for a more ambitious goal, namely to understand the excitation spectrum of molecules on the basis of 
f**«« , the new wave mechanics. They accomplished to exploit the small electron/nucleus mass ratio as an expan- 

■ sion parameter, which then leads to the static Born-Oppenheimer approximation. Since then it has become a 
, standard and widely used tool in quantum chemistry supported by a vast number of mathematical results, the 

V~j ■ classical ones being [6, 12, 18]. 

Beyond excitation spectra and stationary scattering, dynamical processes gain increasingly in interest. Ex- 

■ amples are chemical reactions or the decay of an excited state of a molecule. Such problems are conveniently 
described within a time- dependent version of the Born-Oppenheimer approximation, which is the topic of this 
article. The starting point is the observation that the electronic energy at fixed position of the nuclei serves as 
an effective potential between the nuclei. We call this the zeroth order Born-Oppcnhcimcr approximation. The 
resulting effective Schrodinger equation can be used for both static and dynamic purposes. Of course, the input 
is an electronic structure calculation, which for the purpose of our article we regard as given by other means. 

While there are many physical and chemical properties of molecules explained by the zeroth order Born- 
Oppenheimer approximation, there are cases where higher order corrections are required. Famous examples are 
the dynamical Jahn- Teller effect and the tunneling of singled out nuclear degrees of freedom at a conical inter- 
section of two Born-Oppenheimer energy surfaces. The first order Born-Oppenheimer approximation involves 
geometric phases, which are of great interest also in other domains of non-relativistic quantum mechanics [3]. 
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Our plan is to repeat in some detail the conventional argument for the first order Born-Oppenheimer approxi- 
mation. Only then the reader will appreciate why a more systematic approach is in demand, which will be the 
focus of our contribution. 

Let us start from the molecular Hamiltonian 

ti 2 n 2 

H mo i = -^ A x -- — A y + V(x,y). (1) 

Here x = {xi, . . . , xk} are the positions of the K nuclei and y — {y%, . . . , yjv} the positions of the N electrons. 
The electrons have mass m e and the nuclei have, for notational simplicity, all the same mass m n . V(x, y) is the 
interaction potential, i.e. the Coulomb potential for electrons and nuclei. For the mathematical results to be 
valid in the form we present them, one needs to slightly smear out the charge distribution of the nuclei. This is 
in line with the physical picture that nuclei are not pointlike but extended objects. We transform ([1]) to atomic 
coordinates such that h = 1 = m c and 

*=J^, (2) 

V m n 

which will be our small dimcnsionlcss expansion parameter. Then {l} becomes 

H E = -±e 2 A x +H e (x) (3) 

with the electronic Hamiltonian 

H e (x) = -±A y + V(x,y), (4) 
which depends through V parametrically on x. For complex molecules, one is often forced to model only a 
carefully chosen subset of degrees of freedom. We assume here that the final form is still as in ^ at the expense 
of a suitable modification of H c (x). Implicitly with ((3|) it is required that the initial wave function has a kinetic 
(and hence also a total) energy which is bounded independently of e. Thus the nuclei move slowly. 
The electronic structure problem is the eigenvalue equation 

H e (x)xj[x) = E j (x)x J {x) (5) 

with Xj £ , the Hilbert space for the electronic degrees of freedom. Here "f" is supposed to remind of fast. 
Since electrons are fermions, in our case 7i{ = S a L 2 (M. 3N ) with S a projecting onto the antisymmetric wave 
functions. The eigenvectors in ([5|) are normalized as (Xj( x ) ? Xj'{ x ))Hf = Sjj 1 with respect to the scalar product 
in Tif. Note that the eigenvectors are determined only up to a phase $j(x). Their smooth dependence on x 
will be addressed with more care below. Generically, in addition to the bound states, H e (x) has continuous 
spectrum. 

We label the eigenvalues in ([5]) as 

Ei{x) < E 2 {x) <... , (6) 
including multiplicity. The graph of Ej is the j-th Born-Oppenheimer energy surface. As a rule they have a 
complex structure with many crossings and avoided crossings. Let ip(x) be a nucleonic wave function, ip € Tt s , 
with "s" reminding of slow. For simplicity we take TL S = L 2 (M. 3K ), remembering that to impose the physically 
correct statistics for the nuclei requires extra considerations [23]. States with the property that the electrons 
are precisely in the j-th eigenstate are then of the form 

i>(x)xj(x,y) , (7) 

which we can think of either as a wave function in the total Hilbert space H = Ti. s €3 Tit or as a wave function 
for the electrons depending parametrically on x. In the latter case we abbreviate as ip(x)xj(x) G Tlf for each x. 
With this notation the projection operator onto the states of the form ([7]) is given by 



(8) 
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Since the Xj( x V s are orthonormal, Pj is indeed an orthogonal projection. 

In a molecular collision or in excitations through a laser pulse only a few energy surfaces take part in the 
subsequent dynamics. Thus we take a set / of adjacent Born-Oppenheimer surfaces and call 

P = J2 \Xi(x))<Xi(x)\ (9) 

the projection onto the relevant subspace (or subspace of physical interest). To ensure that other bands are not 
involved, we assume them to have a spectral gap of size a gap > away from the energy surfaces in /, i.e. 

\Ei (x) - Ej (x) | > a gap for all j E I , i E I c . (10) 

Also the continuous spectrum is assumed to be at least a gap away from the relevant energy surfaces. 
If e = 0, then 

e- iHOt Mx)\ Xj (x)) = \xj(x)) . (11) 

However, the Laplacian —e 2 A x weakly couples the x-fibers and P7i is a subspace which is not invariant under 
the true unitary propagator exp[— iH e t\. Still, we take 

PH e P (12) 

as an approximate Hamiltonian for the time evolution in the relevant subspace PTL. To compute PH e P we 
choose (p,ip G PTL, to say ip(x,y) = X^e/ ^j( x )Xj( x jy) with E 7i s and similarly for ip, and sandwich to the 
right and left as 



(<p,Hij;) = J dx<p m (x)*(x m {x), H £ Xn(x))m ipn(x) (13) 

= Y J &%<Pm(x)* 8 mn E m (x) + \e 2 { - S mn A - V • (xm(x), Vx„(x)) Wf 

+ (VXm(x), Xn(x))n t ■ V + (VXm(x), -VXn(a;))w t ) 1/Jn(x) , 

dx4> m {x)*{PH £ P) mn ^ n {x), (14) 

where all derivatives are with respect to x. It is thus natural to introduce the geometric phase, or Berry 
connection, 

A mn (x) = i(Xrn(x) , Vx„(x)) Wf , Amn(x)* = A, lm (x) . (15) 

Here, for each x, Xj( x ),j £ I, is completed to an orthonormal basis through X ji j S P- We also define the 
nuclear momentum operator p — — ieV. Noting that (V^m, Xn) + (x-m, Vx„) = and inserting in the last term 
of (TIB"]) , one obtains 



(PH e P) mn = E m (x)S m n + 5 (P S ^e - eA m l(x)) ■ (pSi n - sAl n (x)) 

eel 

+ \e 2 A m i(x) ■ Ai n (x) . (16) 



PH E P acts on wave functions of the form ip n (x), n E I, i.e. ip E Ti s <8> C^'. 

To make further progress we concentrate on two cases of physical interest. 
1) Let |/| = 1, I = {j}, i.e., we consider a single nondegenerate energy band, which by assumption is isolated 
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from the remaining energy bands. Then Xj{ x ) can chosen to be real and smooth, which implies 

A jj (x) = 0. (17) 

In other words, by a suitable choice of the phase $j(x), see below equation ([5]), i.e., by a suitable gauge, the 
geometric phase can be made to vanish. According to (|13[) this results in 

{PH'P)^ = \ V 2 + E 3 (x) + \e 2 <t>{x) (18) 

with the Born-Huang potential 

4>(x) = <V Xi (z) , (1 - \xj(x))<Xi(x)\) ■ V Xj (x))-H f ■ (19) 

2) Let |/| = 2. For convenience we label asm,n= 0, 1. Eo(x) and E\(x) are allowed to cross. To be very 
specific, and to link to Sectional let us assume that i£l 2 and E o (0) = -Ei(O), while E {x) ^ E\{x) otherwise. 
As first pointed out by Mead and Truhlar [23], if one insists on smoothness of the eigenbasis away from x = 0, 
then xo(x),Xi( x ) are necessarily complex-valued. Amn(x), m,n = 0,1, corresponds to a vector potential with 
a magnetic field which is concentrated at x — and which cannot be gauged away. As to be explained in more 
detail in Section [3l in principle, there are then two choices. In the subspace PH. one can work in the adiabatic 
basis, adopted here, at the expense of having singular coefficients at x — 0. The alternative choice is to pick 
some diabatic basis, which is not an eigenbasis of H e , but has the advantage that in this representation PH £ P 
depends smoothly on x. 

The conventional derivation of the time-dependent Born-Oppenheimer approximation leaves two points in 
the dark. 

(i) Since the total energy is fixed, the nuclei move slowly and thus one has to follow their dynamics over a 
sufficiently long time to see a nontrivial dynamics. In our units the velocities are of order e, hence times have 
to be of order e~ l . At the zeroth order Born-Oppenheimer approximation one neglects terms of order e in 
the intra-band Hamiltonian PH e P, which over times of order e" 1 add up to an error of order 1 in the wave 
function. Thus on this time-scale the first-order Born-Oppenheimer approximation, i.e. the inclusion of the 
Berry connection term, becomes mandatory. But there is a more subtle problem with the standard derivation. 
Since the commutator [P, H £ ] is of order e, it is not clear a priori, cf. (f20|) . if states of the form ([7]) remain even 
approximately of this form for the relevant times of order put differently, if PTL is approximately invariant 
under the time evolution generated by H e for sufficiently long times. Again, this is because naively the order 
e error in the generator may add up to an error of order 1 during times of order e _1 . Thus the substitution of 
H e by PH e P needs justification. 

(ii) If one wants to have the motion of nuclei at higher precision, then the simple ansatz (0) becomes 
questionable. Rather one expects small corrections to the product form. From the way we have presented the 
computation it is not so clear how to include such effects. 

In Section [2] we review mathematical results which settle the two issues raised. The justification of the 
first-order Born-Oppenheimer approximation as in [30] is based on generalizing the standard adiabatic theorem 
of quantum mechanics to a space-adiabatic theorem. The higher order corrections require a systematic scheme 
developed in [5, 8, 22, 24, 26, 28] heavily based on pseudo-differential calculus with operator valued symbols, a 
method which traces back to the pioneering work of Sjostrand [29]. We therefore refrain from explaining the 
general scheme in detail, but instead present the main ideas and a new elementary derivation of the correct second 
order Born-Oppenheimer Hamiltonian. Section [3] deals with applications where the higher order corrections are 
potentially of importance. In particular we discuss the motion near a conical intersection and the reactive 
scattering H 2 + H — > H + H 2 - 
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2. Justification of the Born-Oppenheimer approximation and higher order 

corrections 

In this section we explain several mathematical results concerning the validity of the Born-Oppenheimer 
approximation and higher order corrections to it. We do not restate precise mathematical theorems, which can 
be found in the quoted literature, but focus on structural aspects instead. To this end we first argue how the gap 
in the conventional derivation of the Born-Oppcnhcimcr approximation, as discussed in the introduction, can be 
closed. This understanding will open up the way to systematically determine also the higher order corrections 
to the Born-Oppenheimer approximation. While the general theory describing higher order corrections is rather 
technical, in Section 12.21 we present an elementary computation which yields the second order in s corrections 
to the Born-Oppenheimer approximation. 

2.1. Justification of the time-dependent Born-Oppenheimer approximation 

In the conventional derivation of the zeroth and first order Born-Oppenheimer approximation the effective 
Born-Oppcnhcimcr Hamiltonian is obtained by expanding PH e P in powers of e, where P projects on wave 
functions of the form l[7]). As explained under point (i) above, there is left open an important point in this 
derivation. For PH e P to be a valid effective Hamiltonian in the first place, the dynamics generated by PH e P 
for initial states in P7i must be close to the one generated by the true Hamiltonian H e for sufficiently long 
times, i.e. the difference 

^ e -iH c t/e _ e -iPH c Pt/e^ p 

must be small in an appropriate sense. The Duhamcl formula yields 

r t/e 



(e-^th _ e -iPH°Pt/e^ p = ie -iHH/e J * dg jH's ^pjjep _ jje) ^-iPH^Ps p 

rt/e 

= ie' iHH/e dse iHes {PH £ P - H £ )Pc- lPH ' Ps 
Jo 



rt/e 

= ie- iHH/e dse iHes [P,H 6 ]Pe~ iPH ' Ps . (20) 
Jo 

When acting on wave functions with kintic energy of order 1, i.e. ||eV-0 e || 2 = 0(1), the commutator 



[P,£T]P 



e 2 



P = 0(e) 



contains terms of order e but not smaller, independently of the choice of Xj( x )- Thus in (|20p a naive estimate 
of the right hand side yields an error of order 1. This is because the smallness of the integrand is cancelled by 
the growing domain of integration. Hence it is not clear from a naive perturbation argument that the dynamics 
generated PH £ P and thus the Born-Oppenheimer approximation is close to the true dynamics for sufficiently 
long times. This is directly linked to the question, whether the subspace PH of wave functions of the form 
(J7J) is invariant under the true dynamics for sufficiently long times. We stress this point, because the question 
of invariance of subspaces on which one approximates the dynamics by an effective Hamiltonian is crucial for 
understanding the higher order corrections to the Born-Oppcnhcimcr approximation. 

For the case of the first order Born-Oppenheimer approximation it was shown in [30] how to use arguments 
similar to those used in the proof of the adiabatic theorem by Kato [17] in order to prove approximate invariance 
of PH. In such a space-adiabatic theorem, the positions x of the nuclei take the role played by time t in usual 
iime-adiabatic theory. The mechanism is that the integrand in (|20p is oscillatory and therefore the errors of 
order e in the integrand cannot add up to an error of order 1 even for a domain of integration of size 1/e. The 
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statement proved in [30] (with some typos in the proof corrected in [31]) is the following: Let Pg := l(_oo,B] (H £ ) 
be the spectral projection of H £ on energies below E. Then there exists a constant C < oo independent of e 
such that 



-iH s t/e 



-iPH'Pt/t 



PPf 



B(H) 



<Ce(l + \t\)(l + \E\), 



(21) 



where || • is the norm of bounded operators. 

Thus on subspaces of finite total energy the adiabatic approximation holds with a uniform error of order s. 
Note that (f2"Tj) holds also with Pe replaced by 1(-oo.b] ( — ifA), i- c - f° r finite kinetic energies of the nuclei. For 
large kinetic energies the velocities of the nuclei are no longer small compared to those of the electrons and thus 
the adiabatic approximation breaks down. 

Clearly the first order Born-Oppcnheimer approximation now follows from (|21|) by expanding the adiabatic 
Hamiltonian PH £ P as in (fT3|) . If we define the first order Born-Oppenheimer Hamiltonian as 



H. 



(i) _ 



m.nGl 



\Xm{x))\E m {x)8 mn + \ ^ (p$m£ - eAmi(x)) ■ (pS ln - eA£ n (x))^j (x„{x)\ , 
eei 



where wc omit one term of order e 2 in (11311. then 



-iH e t/e 



PPf 



8(H) 



<Ce(l + \t\){l + \E\ 



(22) 



(23) 



By now it is obvious why we cannot improve the approximation simply by adding the e 2 terms in the expansion 
of PH £ P to the effective Hamiltonian. Their contribution is of the same order as the error in the adiabatic 
approximation (I2ip . One might hope that the order of the error in (I21| can be improved by a more careful 
analysis. This is however not the case, as can be seen directly from the proof in [30] or, alternatively, from 
the scheme to be presented in the following. In particular there is no reason to expect that fp~3|) is the correct 
second order Born-Oppenheimer Hamiltonian and indeed, it isn't. 

2.2. Corrections to the time-dependent Born-Oppenheimer approximation 

But then the question arises whether and in which sense is the first order time-dependent Born-Oppenheimer 
approximation the leading term in a systematic perturbation expansion? This problem was considered several 
times in the literature and the solutions differ not only with respect to the level of mathematical rigor but also 
the formulation of the result itself is not unique. In the following we explain the approach from [26] which 
heavily relies on methods developed by Martinez, Nenciu, and Sordoni [22,24]. The main idea is to replace the 
subspace PTL by a so called almost invariant subspace P £ TL, which is invariant under the full dynamics to higher 
accuracy than PTL. Because of the very technical character of the general construction, we only sketch the main 
ideas and instead give an alternative derivation of the correct second order Born-Oppenheimer Hamiltonian not 
using the pseudo-differential operator machinery. We shortly comment on other approaches at the end of this 
section. 

Since the range of P is invariant only up to an error of order e, the form ([7]) is not adequate for higher 
order approximations. Instead one has to look for slightly tilted subspaces P £ TL which are invariant under the 
dynamics generated by H £ with smaller errors. Indeed one can find an orthogonal projection P £ close to P 
such that for any n € N there is a constant C n < oo independent of e and t satisfying 



iHH/e _ e -iP'H'P't/e\ pe p 



<C n e n \t\. 

B(H) 

It turns out that the projected Hamiltonian P £ H £ P £ can be expanded in powers of e as 



(24) 



P £ H £ P £ =J2 sjH 3 + 0(e n+1 ) =: H £ {n) + 0{e n+1 ) , 

j=0 
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which together with (|24p implies 



P £ P f 



B{H) 



<C n e n \t\ 



(25) 



In the same fashion the projector P e can be expanded in a Taylor series, which for the purpose of the estimates 
in (|24p and (|25|) may be truncated at the appropriate order. Since the zeroth and first order Born-Oppenheimer 
Hamiltonians are just the leading terms in the expansion of Hf n y one could call the n-th order Born- 

Oppenheimer Hamiltonian. But is an effective Hamiltonian which acts on a e-dependent subspace of the 
full Hilbcrt space 7i, while the usefulness of the standard Born-Oppenheimer approximation comes partly from 
the fact that the effective Hamiltonian acts on wave functions depending only on nuclear and possibly a few 
discrete electronic degrees of freedom. Before we explain how to remedy this shortcoming let us briefly explain 
how to construct the projector P e . 

The construction of the projection P e follows a general scheme outlined by Emmerich and Weinstein [8], 
refined in [5] and [24], and finally applied to the Born-Oppenheimer approximation by Martinez and Sordoni [22]. 
The basic idea is to determine the coefficients Pj in the expansion Pf n \ = P + 53j=i £j Pj order by order such 
that Pf n \ is approximately a projection and commutes with H £ up to terms of order s n+ , Writing Po = P it 
clearly holds that Po is a projection and that the commutator 



Pn,--A 



= k(VP )-(-ieV) + -(AP ) 



is of order e when acting on functions of bounded kinetic energy. One now determines Pj inductively by requiring 
that 



(P(n)) 2 -P 



(n) 



0{e 



n+l\ 



(26) 



PL), & 



0(e 



n+l\ 



and assuming the the analogous statement holds already for P? n _iy While the general construction is most 
conveniently done using the theory of e-pseudodifferential operators, let us explicitly determine P\ from the 
above conditions. Because of 



(Po + ePi) 2 - (P + ePi) = e(PoP + PP - P) + 0(e 2 ) 

we must require that 

Pi = PoPi + PiPo + 0(e) 
in order to make the order e term vanish. And the order e term in 



(27) 



[Po + £P lf H E ] = ie(VPo) • (-ieV) + e[P, H c ] + 0{e 2 ) 

vanishes, if Pi satisfies 

[Pi,H c ] = i(VPo) • (ieV) + 0(e) . (28) 
To focus on the simplest case, from now on we assume that Pq(x) projects onto the eigenspace of a single 
eigenvalue Ej(x), i.e. that we are in the situation of item 1) in the introduction. In this special but most 
important case, the unique solution of (j2"B"|) up to order e is obtained by multiplying l|2"8")) from the left and from 
the right by P resp. (1 — Po). The block-diagonal terms are of order e, since Po(VPo)Pq = (1 — Po)(VPo)(l — 
Po) = 0. For the off-diagonal terms one can invert (H e — Ej) on the range of (1 — Po) and finds that 



P 1 = iP (VP )(Pe - Ej)- 1 ^ - P ) • (ieV) + adj. = i(kV) • P (VP )(P C - - Po) + adj. + 0(e) , 
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which also satisfies the requirement (|27[) . The abbreviation 4- adj. means that the adjoint of everything to the 
left is added. Here and in the following we use the fact that the spectral projection Pq(x) and the reduced 
resolvent (H c (x) — Ej(x))~ 1 (l — Pq(x)) are smooth function of x with values in the bounded operators on the 
electronic Hilbert space Hf. This follows from the assumed smoothness of H c (x) and the gap condition. As a 
consequence the commutator of the momentum operator p = — ieV with any such operator yields only a lower 
order term in e, a fact which will be used several times in the following computations. 

In a next step the operator P^ can be turned into a true orthogonal projection by adding a term of order e 2 , 
see [24], which we denote again by Pny Now the question arises, whether the range of P^ is invariant under 
the dynamics generated by H e up to errors of order e 2 . A naive perturbation argument as in (|20|) will only 
yield an error of order e, which, however, could be improved using again the space-adiabatic approach of [30]. 
Alternatively the results of [22,26] show that instead of (|2"Tj) we now have 



iH e t/e 



P^Pf 



8(H) 



<Ce 2 (l + \t\)(l + \E\). 



(29) 



With (|29[) we are in a position to determine an effective Hamiltonian on the range of Pf^-. by expanding 
P^H £ Pfjs in powers of e and keeping terms up to order e 2 . This improves the conventional first order Born- 
Oppcnheimcr approximation by one order. However, because of the momentum operator appearing in Pi, 
the range of Pf^ is not spanned by wave functions of the form ip(x)Xj(x,y) as in Q. Therefore P^H e P^ 
can no longer be seen as an effective Hamiltonian acting on a nucleonic wave function ip{ x ) alone. But this 
reduction in the degrees of freedom of the state space is the crucial feature which makes the first order Born- 
Oppenheimer approximation so useful. In order to retain this feature also for the higher order Born-Oppenheimer 
approximations we thus need to map the range of P^. unitarily to the space L 2 (U. 3K ) of nucleonic wave functions 
in the case of a simple electronic band or to L 2 (M 3K , C' 1 ') in the case of a group of |/| bands. In [26] we construct 
a unitary operator U e : P 6 TC — > L 2 (M. 3K , C' 7 ') and define the effective Born-Oppenheimer Hamiltonian for the 
group of levels in I as 

J?§ = U e P £ H e P e U e * , (30) 

i.e. the full Hamiltonian H £ is projected onto the almost invariant subspace associated with the levels in / and 
then unitarily mapped to the e-independent space of nucleonic wave functions. Thus (|24p becomes 



-iHH/e _ jjs * e -iH^ t/ejje\ pe p 



B(H) 



<C n e n \t\ 



(31) 



where e~ lH ^ ot ^ E is the effective Born-Oppenheimer propagator of the nuclei within the relevant group of bands 
up to any order in e. The expansion of H^ in powers of e yields at the leading orders the usual zeroth and first 
order Born-Oppenheimer Hamiltonian, i.e. the matrix in the brackets in (|22|) . but in addition also the correct 
higher order terms. 

In order to avoid the technicalities of the general construction, let us determine the correct second-order 
Born-Oppenheimer Hamiltonian for the special case of a singe simple level Ej(x) by an elementary calculation. 
For this case we already computed P £ = Pq + eP\ + 0(e 2 ) with Pq(x) — \Xj( x ))(Xj( x )\ an d 



P 1 = -i(-kV) • P Q (VP )(H C - Ej)- 1 ^ - P ) + adj. =: p ■ B + B* ■ p. 

Here and in the following we abbreviate p = —ie\7 x for the momentum operator. Also for II s we make the 
ansatz U e = Uq + eU\ + 0(e 2 ) where clearly 



U (x) = ( Xj (x)\ 



(32) 



is unitary from PqH to L 2 (M. aK ). Note that at this point the construction of U £ is not unique since the choice 
of basis Xj ( x ) enters. 
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Without loss of generality we write U\ = UqA for some operator A on Ji. Then the requirement that U e * is 
unitary up to 0{e 2 ) yields 



(33) 



(U + eUx){U; + eU*) = 1 + e{U A*U* + U Q AU*) + 0(e 2 ) = 1 + 0{e 2 ) . 
The requirement that the range of U e * is the range of P e up to 0(e 2 ) yields 

(1 - P - ePi){U; + eUl) = e(A* - P A* - P 1 )U^ + 0(e 2 ) I 0(e 2 ) , 

i.e. that 

(1 - P a )A* = PP + 0{e) . 
A solution of this equation is A* — B* • p, which also makes the 0(e) term in (f33|) vanish. Thus 

Ui = U p ■ B 

yields an almost unitary = Uq + sUi that almost intertwines P^H and L 2 (M. 3K ). It is important for the 
following that can be modified by a term e 2 /7| of order e 2 (i.e. £/| is 0(1)) to make it a true unitary exactly 
intertwining and L 2 (M. 3K ). As in the case of P^ we denote the modified operator Uf^ = Uq+sUi +e 2 C/| 

again by the same symbol. It is shown in [26] how to construct f/f explicitly, but its exact form is not 
important for the following. We will only use that unitarity of the modified Uf^ together with the fact that 
IWo + = U PiU$ = implies that 



U x Ut + IWf * + U$US = 0{e) . 

We thus have that 

and therefore with ([29)1 that (|35[) approximates the true time evolution up to errors of order e 2 , 



[e-^t/e _ V e* ^U^P^H'P^U^t/e^ pe^ P E j < C E 2 (1 + \t\) 



(34) 



(35) 



(36) 



Hence we can now expand U^P^H e P^U^ in powers of e to obtain the second order Born-Oppcnheimer 

Hamiltonian, which now acts on the e- independent space L 2 (W? K ) of nucleonic wave functions. The expansion 
yields 

U^P^H'P^U^ = UfaH'Uft = (U + eXJ 1 + e 2 U e 2 ) + P c ) (U* + eU* + e 2 [/f*) + 0(e 3 ) 

= C/o( 2 j+Pc)f/o (37) 
+ e {u Q (£+H^Ut + adj.) (38) 

+ £ 2 [/l(^+Pe)t/l* 



e 2 (U ^+P e U. 



(39) 
(40) 
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We evaluate the four terms l|5Tj) - (|4"0")) separately. Expanding (|57|) yields 

U + fle) US = I T + E ] +U Q [^-,U* 

= ^ + Ej - ieU {VU*) ■ p - £U (AU£) 



i(p-e^) 2 + S J + 40, (41) 



where we abbreviated 



A(x)=is( Xj (x),V X j(x)) (42) 



for the Berry connection coefficient and 

cf>(x) = (V Xj (x),-(l-P )V Xj (x)) (43) 

for the Born-Huang potential. This is, as expected, exactly the usual Born-Oppenheimer Hamiltonian (|13D . 
However, as explained before, in order to obtain the correct second order Born-Oppenheimer Hamiltonian, we 
also must take into account the terms (|55|) -([4T) ]) which stem from the fact, that not PTL but is the correct 

adiabatically invariant subspace at that order. Expanding ([55} yields 



o 



Uo(£+H^UZ +adj. = U ^B* -p + p- B^U* + U {H C B* -p + p- BH C )U { 

= U ([£,B*].p + p-[B,£])uS 

= -isU Q (p- VB* -p-p-VB -p)Ua + 0{e 2 ) 

= -2eU oP - P (VP Q )(H C - E J )- 1 (\7P Q )P ■ P U^ + 0(e 2 ) 

=: -2eM(x,p) + 0(e 2 ), (44) 

where we used that VBP = iP (\7P )(H c - E J )- 1 {\7P )P and P VB* = -iP (VP )(H e - E j )- 1 (VP )P . 
Recall that the differences in operator ordering between p = — ie\7 x and operators depending on x are of lower 
order in e. So any other operator ordering for the quantization of A4(x,p) than the one used in (|38|) works as 
well at the given order. Put differently, any quantization rule for the symbol (=function) 

M(x,p) = U oP -Po(VPo)(H c -E j y 1 (VPo)P -pUS = <V Xj -p, (H e - Ej)- 1 ^ - P Q )p ■ V Xi > Wf (45) 

will do the job. The simplest symmetric choice for M. is presumably 

3^ 

(Mip)(x) = ^2 2 ( m ek(x){-ied Xl ,){-ied Xk ) + (-ied Xe )(-isd Xk )m ek (x)jip(x) , 
e,k=i 

where m is the x-dependent matrix 

m ek (x) = (dtXj{x),{H c (x) - Ej(x))~ x (l - P {x)) d k Xj(x)) Hf 

For {33]) we find 

Ux{£+H e )U{ = U p-B{£ + H D )B* - P U* a 

= U p- B{H C - E,)B* -pU* + U oP - B{£ + E 6 )B* ■ p U* 

= U oP - Po(VP )(#c - E j )- 1 (VP )P -pUZ + Uop- BB* ■ pU£ + Ej) + 0(e) 

= M(x,p) + U 1 U?{£+E j ) + G(e). (46) 
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Finally the term (|4"0|) cancels the second remaining term in : 

U + Be) U! * + U! + H c ) [7 * = 17 (£ + Ej) Ul * + C/f + E 3 ) [/* 

= (U Ul * + UZUS)(£ + E j ) + 0(e) 

= -UxUUi+E^ + Oie), (47) 

where we used ([34|) . Collecting all the results we find for the second order Born-Oppcnhcimcr Hamiltonian of 
a simple isolated eigenvalue band Ej that 

#bo = W - £*4(z)) 2 + Ej(x) + £<f>(x) - e 2 M{x lP ) , (48) 



where we recall that p = — \e\I x and that A(x), <f)(x) and A4(x,p) are defined in ([42]) , (|43|) and (|45|) . 

The relation between the dynamics generated by the Born-Oppenheimer Hamiltonian -fffj an d the true time 
evolution is now 

(e-WH/s _ ue* e -iH| t/ S[/ eJ pe ^ p E ^ < Ce 2 (1 + |^l) ■ (49) 

Put differently, one can construct approximate solutions v E , bo(^) to the full molecular Schrodinger equation 

ied t y(t,x,y) = H £ ^{t,x,y), 6 L 2 (R 3K ) <g> W f (50) 

by solving the effective Born-Oppenheimer Schrodinger equation for the nuclei only 

ied t m = fliolK*) , 6 i 2 (K 3A ') , (51) 

and defining "FboW = U £ -^ip(t). The difference between \I/Bo(i) and the true solution of (fSTj) with the same 
initial condition ^(0) = 'I'bo(O) G PJ^PeH is of order e 2 in the norm of K. However, there are many cases 
where it suffices to know tp(t) and it is not necessary to map it back to a full molecular wave function \&bo(£)- 
For example for the position distribution of the nuclei one has that 

(*Bo(M),*Bo(M))« f = (U e {1 p(t,x),U e {1 ]ip(t,x)) n{ = \^(t,x)\ 2 + 0(e 2 ). 

We continue with several remarks on the second order Born-Oppenheimer Hamiltonian. The correct second 
order Born-Oppenheimer Hamiltonian (|48[) contains the additional term e 2 A4(x 1 p) as compared to the con- 
ventional expression containing terms of that order, see (fT8"|) . (|19[) . In contrast to the Born- Huang potential 
(f>(x), M(x,p) introduces a velocity-dependent correction in the form of an ir-dcpcndent effective mass tensor 
m. In Section [3] we compute i? BO explicitly near a conical intersection of electronic energy levels, where the 
corrections to the conventional first order Born-Oppcnhcimcr approximation become important. A further issue 
is the gauge invariance of -ff BO - The Berry connection coefficient A(x) clearly depends on the choice of Xj( x ) 
and, in our particular context, can be made to vanish by a suitable choice of phase. However, although defined 
through Xj( x )i the terms of order e 2 are gauge invariant, as can be seen from writing them only by means of 
the projection Pq(x), 

<t>(x) = Tr nt (yP (x)- VP (x) (1 - P (x)) , (52) 
M(x,p) = TT H{ ^(p-VPo(x)) 2 (H c (x)~E ] (x)y 1 (l-P (x))^ . (53) 

In general, the structure of i? BO is that of a semiclassical Hamiltonian for the nuclei, since the momentum 
operator still carries the small parameter s. Therefore it is natural to study in a second step the semiclassical 
limit of the nuclear Schrodinger equation (fSTj) . There is a variety of methods available for such a semiclassical 
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analysis, most prominently the WKB approximation, scmiclassical wave packets and Wigner functions. All 
these techniques have been applied to further study the standard Born-Oppcnhcimcr approximation and can 
also be used to investigate the higher order Born-Oppenheimer approximations. But since the semiclassical 
limit is conceptually and mathematically different from the adiabatic approximation considered here, we do not 
further comment on it. 

We conclude this section with a few remarks on the literature. In the physics literature the correct second 
order Born-Oppenheimer Hamiltonian (|48|) was first obtained by Weigert and Littlejohn in [32] for the case 
of matrix valued H c (x). They approximately diagonalize the Hamiltonian H E on the full space Tt, while we 
first reduce to an appropriate adiabatic subspace corresponding to the electronic levels of interest and then 
approximate the Hamiltonian on that subspace. This has the advantage that we can allow for H c (x) having 
continuous spectrum or level crossings outside the spectral part of interest. 

Another mathematical approach to the time-dependent Born-Oppenheimer approximation, historically the 
first one, is due to Hagedorn and later Hagedorn and Joye [11,13,14,16]. In their approach the goal is not 
to construct an approximate Born-Oppenheimer Hamiltonian but to find directly approximate solutions of the 
molecular Schrodinger equation. This is achieved by approximating the nucleonic wave function by localized 
semiclassical wave packets. Hence, in the Hagedorn/Joye approach the adiabatic and the semiclassical approx- 
imation are done in one package. 



The aim of this section is twofold. In the first part we show how the previous scheme extends to a family 
of energy bands, yielding a multiband effective Hamiltonian. While this effective Hamiltonian is the starting 
point for an analysis of the propagation of the wavefunction through eigenvalue crossings, we will not address 
this interesting problem here and instead refer the reader to the literature [9,10,15,19], see also [7] for a related 
result in a time-independent setup. 

In the second part, Section 13. 2\ we return to the one-band Born-Oppenheimer dynamics and study in a 
simple case the behavior of the Berry connection A and of the second-order corrections tfi and M. near a conical 
crossing. We argue that this behavior is, in a sense, universal and leads potentially to observable effects in 
chemical exchange reactions. 

3.1. The multiband effective Hamiltonian 

We consider a family of t eigenvalue bands {Ej}j e i, I = £, which may cross each other but which are 
separated by a gap from the rest of the spectrum, see ([TO]) . We denote as Po(x) the eigenprojector of H c (x) 
corresponding to such a family of bands, with dimRanPo^) = 

By the construction outlined in the previous section, to such a family of bands there corresponds an almost- 
invariant subspace. More precisely one constructs an orthogonal projector P £ = Pq + 0(e) satisfying (|24[). 
In order to describe the dynamics inside the almost-invariant subspace in a physically transparent way, one 
constructs a unitary operator U 6 which intertwines RanP £ and a reference space, i.e. 



3. Dynamics near conical crossings 



U e : RanP £ 



Href : 



L 2 (R 3K ,C e )- 



(54) 



To the lowest order U e = U Q + 0(e), with 



t 




(55) 



a=l 



where {e a } a =i,...,e is the canonical basis in C m and {<p a (x)}a=i,...,£ is any orthonormal basis spanning RanPo(a;) 
and depending smoothly on x. The freedom in choosing such a basis is an additional feature with respect to the 
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previous section. If \jPa}a=x,...,l is anv such basis, by the construction in [26] one obtains a self-adjoint operator 
i?B , such that (|3 1 [) is satisfied. It has the asymptotic expansion 



J?io := U*P°H°P°U° 

\p 2 + W(q) + f (p • A(x) + A(x) ■ p) + 0(e 2 ), 



where we introduced the £ x £ matrices 

W{x)ab = (<Pa(x), H e {x)<Pb(x)) Ut , 

(57) 

A{x) ab = i((p a (x),Vtpb(x)} Hi , a,b=l,...,£. 

The eigenvalues of W(x) are the energies {Ej(x)} j<zi of the electronic Hamiltonian H c (x). The expansion in 
the second line of (|56[) is the generalization of (|2"2"|) to the case of an arbitrary basis. 

Since H Cl as defined in ([4J, satisfies time-reversal symmetry the eigenfunctions {<p a }a=i,...,£ can be chosen 
to be real valued. With this restriction W(x) becomes a real symmetric matrix and iA{x) becomes a real 
antisymmetric matrix, in particular with vanishing diagonal elements. 

The effective Hamiltonian (|56p is clearly not unique. Indeed, one can consider a different basis {<f a }a=i....,e, 
with 

i 

lfia(x) = ^ G{x)ab ftb(x) 
6=1 

where G(x) is some orthogonal matrix depending smoothly on x. The elements of the effective Hamiltonian in 
these two basis are related by 

W{x) = G{x)- 1 W(x)G(x), (58) 
A(x) = G(x)- 1 A(x)G(x) + G(x)~ 1 VG(x). (59) 

Notice that, even with time-reversal symmetry, in general the Berry connection A cannot be removed by a 
change of basis (a non-abelian change of gauge). Indeed, by (|59[) this is possible if and only if there exist a 
smooth solution G(x) of the system of differential equations 

diG(x) = -Ai(x) G{x), i = l,...,3K, 

which as a necessary condition requires the vanishing of the corresponding matrix valued field 

u^ij — i (d^Aj djAi) ~\~ AjA[ AiAj. 



The non-uniqueness of the effective multiband Hamiltonian (|56p is irrelevant as far as expectation values of 
physical observables are concerned, since a change of basis affects both the Hamiltonian operator and the 
operator representing the observable quantity. However, the dependence on the basis may be relevant if one 
considers H^ as a bona fide model Hamiltonian to investigate the dynamics of the wavefunction in a family of 
bands. We will comment on this point in the following. 

3.2. Higher-order corrections near conical crossings 

We come back to the one-band Born-Oppenheimer dynamics, i.e. I = {j}, but consider now the case when 
the relevant energy band Ej crosses other bands, so that (ITU1) does not hold. The dynamics of a wavepacket 
localized far away from the crossing points, and in the almost-invariant subspace corresponding to the relevant 
energy band Ej, is approximately governed by the Born-Oppenheimer Hamiltonian (|48[) . (While in Section [2] 
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we discussed only the case of a globally isolated energy band, the more general case has been studied in 
[26,31]). Since the adiabatic approximation breaks down near conical crossings one expects that the higher 
order corrections A, </> and M. diverge near such crossings. It has been pointed out by Berry and Lim [2] that 
the Born-Huang potential <f> creates a diverging repulsive force at the crossing points, and some authors argued 
that the wavefunction must vanish at those points. While it is certainly true that <f> leads to a repulsive force, 
we emphasize that a complete analysis should take into account all contributions coming from the second order 
corrections. 

Since a crossing involves generically just two bands, we focus on a family of £ — 2 bands, E±, separated by a 
gap from the rest of the spectrum. We also focus on the case of conical crossings [15], which are, in a sense, the 
simplest non-trivial kind of crossings. The set of conical crossing points M con C IR 3 ^ is generically a manifold 
of codimension 2, see [25]. Since we are interested in the dynamics in the transverse direction we immediately 
assume x G M 2 , x = being the only conical crossing point. The conical structure at crossing means that 
E±(x) = ±C±\x\ as |z| -> 0. 

The dynamics for the family {E+,E-} is approximately described by an effective Hamiltonian in the form 
(|56p . The choice of the basis appearing in (|55j) requires particular care near conical crossing points. As pointed 
out already, away from the crossing manifold it is convenient to use the adiabatic basis, i.e. a basis {xa}a=± 
consisting of eigenfunctions of the electronic Hamiltonian, see ([S]). Such a basis is uniquely defined, up to a 
choice of the phases, away from the crossing manifold. Clearly, with this choice, W(x) becomes a diagonal 
matrix. On the other hand, Vx(x) generically diverges as x approaches a conical crossing, and so does the 
Berry connection term, yielding an effective Hamiltonian with singular coefficients. Diverging quantities are 
unstable when numerical discretization schemes are employed. Thus, if one is mainly interested in the dynamics 
near the crossing points, it is advantageus to skip condition ([5]), and to work in a different basis, with the 
property that it depends smoothly on x. Such a basis is called diabatic. 

The next step is a reasonable truncation of the e-expansion appearing in (|56|) . We are interested in the one- 
band Born-Oppenheimer dynamics near a conical crossing, where the functions A+, <fi+ and M + are singular. 
Thus in (|56p we retain the terms of order 0(1) and between the terms of 0(e) or higher we retain only the 
singular contributions. Therefore, if (|56p is written in a diabatic basis, so that A{x) is a smooth function, we 
truncate the expansion in (|56|) at the leading order, obtaining the model Hamiltonian i/dia = \p 2 + W(x) (in 
a suitable diabatic basis). Viceversa, if the adiabatic basis is used, the singular term p ■ A + A ■ p has to be 
included. 

The singular behavior of A+ at the conical crossing does not depend on the choice of the basis in (f55|) . 
Indeed, in the adiabatic basis W(x) is diagonal, and the one-band Berry connection for the upper band is given 
by A&dm.+ = i (x+i ^X+)t-c { - On ^ ne other side, in a generic basis {</? Q }a=i,2 one has a non-diagonal W(x) with 
eigenvectors £±(2;) £ C 2 . In this basis the 0(e) term contains 

Adi^+ix) = i(t+{x), A(x)Z + (x)) c2 +i(Z + (x),VZ+(x)} c2 , (60) 

where A(x) a b = i (</? (x), V x iy5b(x)) Wf . The first term in (|60p corresponds to the reduction from H e to the 
2-band Hamiltonian (|56p . while the second term is due to the reduction to a specific energy band, namely E + . 
Note that in the diabatic representation the singular contributions come from the eigenvectors of W , while 
the first term in (f6"0"| is smooth and can thus be neglected for our purposes. The two basis are related by 
Xa(x) — G a b(x)(pb{x), with G(x) an orthogonal matrix smoothly depending on x outside the crossing point. 
With the help of (|5"5)) and (fST))) one easily checks that Alia, + = Aidia, + • A similar behavior is expected for the 
second order terms 4> + and M+, although we are not in position to discuss this point here. 
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With the purpose of having an explicit example, we now focus on the model Hamiltonian -f/dia = — \p 2j rW(x), 
acting in L 2 (R 2 ,C 2 ), with x — (xi,x 2 ) — (\x\ costp, \x\ sintp) and 

W(x)=c( Xl X2 )=C\x\( COSip S ' mip V 0. (61) 
\ X2 —x\ J V sm(p —coscp J 

We freely switch between polar and cartesian coordinates in the following. The eigenvalues of W(x) are E±(x) = 
±C |x| and a smooth family of eigenfunctions for x ^ is given by 

(*) = e-/ 2 ( I s ! ) £_(x) = e-/ 2 ( -Jl n # ) . (62) 

A direct computation yields 



sm -k J \ cos -k 



Vi-(.r) ^-(-^ + (x)+i^(x))e lp , (63) 



where e v = |x| {x2, —x\) in cartesian coordinates. The derivatives V£± are divergent at the crossing point. 
One easily computes the higher order Born Oppenheimer approximation for the upper band, with the result 



A+(x)=i^ + {x),V x ^+{x)) = -J^e v , (64) 



Using 



one obtains 



<t>+{x) = <V^+, •(! - P+(x))V x £+) = (65) 

&\x\* 

(W(x) E+ix))- 1 (1 - P + {x)) = (1 - P + (x)) 

1 2 1 

M+ = ^r- Pi(d^ + (x),(l-P + (x))d^ + (x))p j = L(x, P )—- E L(x,p), (66) 

1 1 i,3'=l 1 

where L(x,p) = x\p 2 — £2Pi is the angular momentum operator. Thus the effective Born-Oppenheimer Hamil- 
tonian to second order reads 

K = \Cp+ ^) + C \x\ + £ 2 -L + e*L(x,p)^L(x, P ) + 0(e% 



For the Hamiltonian of the lower band only A4- = —A4+ changes, while A- = A+ and 0_ = <p+. The inverse- 
square potential corresponds to the repulsive force of Berry and Lim. On the other side, M. contributes to the 
effective Hamiltonian with a more singular term, whose sign depends on the electronic state: on the upper band 
it leads to a repulsive force, while on the lower band it is of the same order, but attractive. In both cases it 
may dominate the contribution coming from the Born-Huang potential. 

It is instructive to reexpress -ffdia in the basis (f6"2"| . Denoting by S the unitary operator corresponding to the 
change of basis, one gets 

nrr o-l 1 2 ( E+(x) \ £ / 1 i \ S 2 { 1 i \ 
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Since the transformation is unitary, no information as compared to -ffdia has been discarded. If for the dynamics 
in the upper band one merely restricts -ff a di to a, b = + , one recovers the Born- Huang potential but misses the 
A4 correction. 

So far we focused on the specific model (|6ip . One may ask wether the singularities computed in (|6"4")l . (|b3|) 
and (|66p are independent from the specific Hamiltonian considered. As for the Berry connection A+, we already 
pointed out that the singularity does not depend on the choice of the diabatic basis used to compute W{x). 
Moreover, the sign of M. does not depend on the specific model. Indeed if two bands come very close, with 
E n (x) < E n +\{x), but themselves are well separated from the others bands, the main contribution to the 
reduced resolvent comes from the n-th band. Thus 

-(W(x)-E n+1 (x)y 1 P r t +1 (x) = - {E m {x)-E n+l {x))- l Pi{x) 

» -(E n (x)-E n+1 (x)y 1 P^(x)>0. 
Thus yVfn+i is positive and, by the mirror argument, A4 n is negative. 

3.3. Second order corrections in chemical reactions 

The search for observable effects of the Berry phase has been a wide field of investigations [3] . In the same 
spirit, one may ask if there are observable effects related to the second order corrections 4> and M. A good 
candidate is the chemical exchange reaction involving a system of three hydrogen atoms: H 2 + H — > H + H 2 - 
Indeed, the intermediate transition molecule H3 exhibits an interesting structure: two electronic energy bands 
form a family of conical intersections, corresponding to those highly-symmetric configurations where the nuclei 
are at the vertices of an equilateral triangle. 

Experimental studies of H2 + H (which are usually performed by using one of its isotopic analogues, as 
D + H-2 — > DH + H, so that the reactants are labelled) have reached a very high level of precision. For example, 
in [1] the rates are measured for reactions in which both the initial and the final hydrogen molecules are in a 
specific vibrational and rotational state. In fact, one expects the state-to-state cross section to be more sensitive 
than the total cross section to phase interference of the Born-Oppenheimer wave function. 

The basic interference mechanism is simple. The labelled reaction A + BC — > AB + C may happen in two 
qualitatively different ways: either (i) the incoming atom A binds directly with B, or (ii) atom A first approaches 
atom C before finally binding with atom B. In the Born-Oppenheimer approximation, this two possibilities 
correspond to classical paths on opposite sides of the conical crossing. Thus quantum interference is expected. 

On the theoretical side, Kupperman and Wu [33] computed state-to-state differential cross section for the 
H + H% reaction, showing significant difference with previous calculations, where the Berry phase was neglected. 
Comparison between theory and experiments was firstly done for D + H2 — > DH + H. In early computations, the 
inclusion of the Berry phase improved significantly the agreement between theory and experiments. However, 
it has been later pointed out, the agreement could be due to accidental cancellations of errors, one source 
of error being the inaccuracy of the electronic structure data in the high-energy region. The results arc still 
controversial, and no general consensus has been reached [3]. 

While we do not enter in the controversy, we comment on the relevance of second-order Born-Oppenheimer 
corrections in this kind of reactions. First of all, some orders of magnitude: in the identical atoms system the 
conical intersection exists at E cross s» 2, 7 eV above the minimal energy of the H2 molecule, and the minimal 
potential barrier is Eq 0.43 eV, corresponding to the collincar H—H—H configuration. We focus on a situation 
in which the kinetic energy is smaller than E cross , so that diabatic transitions to the upper band are negligible, 
but sufficiently larger that Eq, so that the probability that the wavefunction follows a path of kind (ii) is not 
negligible. In this energy interval the wavefunction is likely to explore regions of the nucleonic configuration 
space close to the crossing point. The second order corrections <f> and A4 are singular at the crossing point, 
and thus certainly relevant in the nearby region. Although a careful analysis requires quantitative estimates, 
we conclude that the second-order terms might be relevant in this kind of reactions. 
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